Simulations fitness data
Laure Olazcuaga, Nicolas Rode
04-February-2021
1 Colors plot
colors <- c("LMER" = "#DCA237", "Anova" = "#459B76")
lines <- c("Balanced" = "solid", "Unbalanced" = "dotted")2 BALANCED DESIGN
3 Normal distribution
3.1 Simulation
# Power Analyses
nb_simul <- 10000
#Perform siluations
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "normal",
rho = -0.05, rho_ng = -0.05,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
#Add indic SA
sim <- add_indic_sign(sim)
#Add number of simul for each SA value to keep only 100 simuls per SA value
sim <- add_sim_number_SA(sim)
## Local adaptation (genetic)
tapply(sim$IndicGen[sim$index_SA_Gen<101], sim$SA_Gen[sim$index_SA_Gen<101], length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 100 100 100 100 100 100 100 100 100 100 100 37 19 1tapply(sim$IndicGen_aov[sim$index_SA_Gen<101], sim$SA_Gen[sim$index_SA_Gen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.6000000 0.5800000 0.6200000 0.5200000 0.6600000 0.6400000 0.6700000 0.7000000
## 0.9 1 1.1 1.2 1.3 1.4 1.5
## 0.6400000 0.7800000 0.7600000 0.8300000 0.8648649 0.8421053 1.0000000## Local adaptation (non genetic)
tapply(sim$IndicNonGen[sim$index_SA_NonGen<101], sim$SA_NonGen[sim$index_SA_NonGen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.0000000 0.1100000 0.1800000 0.2000000 0.2300000 0.3100000 0.2800000 0.5400000
## 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 0.4900000 0.5400000 0.6600000 0.6938776 0.6346154 0.9230769 0.8000000 1.0000000tapply(sim$IndicNonGen_aov[sim$index_SA_NonGen<101], sim$SA_NonGen[sim$index_SA_NonGen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.0000000 0.1400000 0.1800000 0.1900000 0.2300000 0.2900000 0.2700000 0.5400000
## 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 0.4900000 0.5300000 0.6600000 0.7040816 0.6346154 0.9230769 1.0000000 1.0000000tapply(sim$IndicGen, sim$SA_Gen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 107 496 1127 1633 1919 1672 1307 868 464 234 106 37 19 1tapply(sim$IndicNonGen, sim$SA_NonGen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 4 111 493 1054 1759 1912 1695 1232 853 451 267 98 52 13 5 13.2 Plot
####################################
###### Plot genetic
####################################
#Number of simuls for each SA value
tapply(sim$index_SA_Gen, sim$SA_Gen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 107 496 1127 1633 1919 1672 1307 868 464 234 106 37 19 1#Save the same number of simuls for each SA value
nb_per_SA=100
sim_select_genetic <- select_sample_simul(dataset_sim = sim,
type_SA = "genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#To check
# dim(sim)
# dim(sim_select_genetic)
# tapply(sim_select_genetic$index_SA_Gen, sim_select_genetic$SA_Gen, length)
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_genet = sort(unique(sim_select_genetic$SA_Gen)),
prop_val_sign = tapply(sim_select_genetic$IndicGen,
sim_select_genetic$SA_Gen, mean),
prop_val_sign_aov = tapply(sim_select_genetic$IndicGen_aov,
sim_select_genetic$SA_Gen, mean))
#Plot
plot_genetic_normal<-ggplot2::ggplot(data = data_prop, aes(x=val_SA_genet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 4) +
geom_line(aes(y=prop_val_sign, color="LMER", linetype = "Balanced"), size = 1.5) +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova", linetype = "Balanced"), size = 0.75) +
ylab("Proportion of simulations with\nsignificant genetic local adaption") +
xlab("Intensity of genetic local adaptation (SA)") +
ggtitle("Normal data") +
theme_LO_sober +
ggtitle("Normal data") +
labs(color = "Methods", linetype = "Design") +
scale_color_manual(values = colors) +
scale_linetype_manual(values = lines)
plot_genetic_normal
# #Save plot
# name_plot<-paste0("Simul_powertest_genetic_normal_",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# plot_genetic_normal, base_height = 10/cm(1), base_width = 14/cm(1), dpi = 1200)
#################################
### Plot non-genetic
#################################
#Number of simuls for each SA non genetic value
tapply(sim$index_SA_NonGen, sim$SA_NonGen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 4 111 493 1054 1759 1912 1695 1232 853 451 267 98 52 13 5 1#Save the same number of simuls for each SA value
nb_per_SA ## [1] 100sim_select_nongenetic <- select_sample_simul(dataset_sim = sim,
type_SA = "non-genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#To check
# dim(sim)
# dim(sim_select_nongenetic)
# tapply(sim_select_nongenetic$index_SA_Gen, sim_select_nongenetic$SA_Gen, length)
# tapply(sim_select_nongenetic$index_SA_NonGen, sim_select_nongenetic$SA_NonGen, length)
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_nongenet = sort(unique(sim_select_nongenetic$SA_NonGen)),
prop_val_sign = tapply(sim_select_nongenetic$IndicNonGen,
sim_select_nongenetic$SA_NonGen, mean),
prop_val_sign_aov = tapply(sim_select_nongenetic$IndicNonGen_aov,
sim_select_nongenetic$SA_NonGen, mean))
colors<-c("LMER" = "#DCA237", "Anova" = "#459B76")
#Plot
plot_nongenetic_normal<-ggplot2::ggplot(data = data_prop, aes(x=val_SA_nongenet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 4) +
geom_line(aes(y=prop_val_sign, color="LMER", linetype = "Balanced"), size = 1.5) +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova", linetype = "Balanced"), size = 0.75) +
ylab("Proportion of simulations with\nsignificant non-genetic local adaption") +
xlab("Intensity of non-genetic local adaptation (SA)") +
theme_LO_sober +
ggtitle("Normal data") +
labs(color = "Methods", linetype = "Design") +
scale_color_manual(values = colors) +
scale_linetype_manual(values = lines)
plot_nongenetic_normal
#Save plot
# name_plot<-paste0("Simul_powertest_nongenetic_normal_",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# plot_nongenetic_normal, base_height = 10/cm(1), base_width = 14/cm(1), dpi = 1200)
#################################
### Common plot
#################################
power_plot_normal<-cowplot::plot_grid(plot_nongenetic_normal,plot_genetic_normal,
labels=c("A", "B"),
label_size = 15,
ncol =1, nrow = 2,
# hjust = 0, vjust = 1,
scale = c(1, 1))
power_plot_normal
#Save plot
name_plot<-paste0("Simul_powertest_normal_",
nb_simul,"nbsimul_",nb_per_SA,
"nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# power_plot_normal, base_height = 20/cm(1), base_width = 14/cm(1), dpi = 1200)3.3 False positive rate analyses
### No genetic and no non genetic effects
nb_simul = 200
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "normal", rho=0, rho_ng=0,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
mean(sim$IndicGen)## [1] 0.01mean(sim$IndicGen_aov)## [1] 0.01## Local adaptation (non genetic)
mean(sim$IndicNonGen)## [1] 0.02mean(sim$IndicNonGen_aov)## [1] 0.025### Genetic effects and no non-genetic effects
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "normal", rho=-0.1, rho_ng=0,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
tapply(sim$IndicGen, sim$SA_Gen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.00000000 0.16666667 0.06666667 0.35000000 0.15384615 0.13333333 0.30769231
## 0.9 1 1.1 1.2 1.3
## 0.26086957 0.54545455 0.66666667 1.00000000 1.00000000tapply(sim$IndicGen_aov, sim$SA_Gen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.00000000 0.16666667 0.06666667 0.35000000 0.15384615 0.13333333 0.30769231
## 0.9 1 1.1 1.2 1.3
## 0.26086957 0.54545455 0.66666667 1.00000000 1.00000000## Local adaptation (non genetic)
mean(sim$IndicNonGen)## [1] 0.02mean(sim$IndicNonGen_aov)## [1] 0.025### Non genetic effects and no genetic effects
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "normal", rho=0, rho_ng=-0.1,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
mean(sim$IndicGen)## [1] 0.08mean(sim$IndicGen_aov)## [1] 0.08## Local adaptation (non genetic)
tapply(sim$IndicNonGen, sim$SA_NonGen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## 0.0000000 0.1250000 0.2000000 0.2187500 0.4285714 0.2702703 0.5000000 0.5333333
## 1 1.1 1.2
## 0.2500000 0.5000000 1.0000000tapply(sim$IndicNonGen_aov, sim$SA_NonGen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## 0.0000000 0.2500000 0.1500000 0.2187500 0.4285714 0.2702703 0.5000000 0.5333333
## 1 1.1 1.2
## 0.2500000 0.5000000 1.00000004 Poisson distribution
4.1 Simulation
nb_simul <- 10000
#Perform siluations
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "poisson", rho = -0.05, rho_ng = -0.05,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
#Add indic SA
sim <- add_indic_sign(sim)
#Add number of simul for each SA value to keep only 100 simuls per SA value
sim <- add_sim_number_SA(sim)
## Local adaptation (genetic)
tapply(sim$IndicGen[sim$index_SA_Gen<101], sim$SA_Gen[sim$index_SA_Gen<101], length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 100 100 100 100 100 100 100 100 100 100 100 37 19 1tapply(sim$IndicGen_aov[sim$index_SA_Gen<101], sim$SA_Gen[sim$index_SA_Gen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.6000000 0.6400000 0.6600000 0.6100000 0.7200000 0.7000000 0.8300000 0.8400000
## 0.9 1 1.1 1.2 1.3 1.4 1.5
## 0.8300000 0.8900000 0.8900000 0.9400000 0.9189189 0.9473684 1.0000000## Local adaptation (non genetic)
tapply(sim$IndicNonGen[sim$index_SA_NonGen<101], sim$SA_NonGen[sim$index_SA_NonGen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.0000000 0.0200000 0.1100000 0.1700000 0.3000000 0.4300000 0.4300000 0.6000000
## 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 0.6100000 0.6800000 0.8100000 0.8469388 0.7884615 0.8461538 0.8000000 1.0000000tapply(sim$IndicNonGen_aov[sim$index_SA_NonGen<101], sim$SA_NonGen[sim$index_SA_NonGen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.2500000 0.1300000 0.2000000 0.2700000 0.3600000 0.4200000 0.4300000 0.5700000
## 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 0.6400000 0.6800000 0.8200000 0.8673469 0.7884615 0.8461538 0.8000000 1.0000000tapply(sim$IndicGen, sim$SA_Gen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 107 496 1127 1633 1919 1672 1307 868 464 234 106 37 19 1tapply(sim$IndicNonGen, sim$SA_NonGen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 4 111 493 1054 1759 1912 1695 1232 853 451 267 98 52 13 5 14.2 Plot
####################################
###### Plot genetic
####################################
#Number of simuls for each SA value
tapply(sim$index_SA_Gen, sim$SA_Gen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 107 496 1127 1633 1919 1672 1307 868 464 234 106 37 19 1#Save the same number of simuls for each SA value
nb_per_SA=100
sim_select_genetic <- select_sample_simul(dataset_sim = sim,
type_SA = "genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#To check
# dim(sim)
# dim(sim_select_genetic)
# tapply(sim_select_genetic$index_SA_Gen, sim_select_genetic$SA_Gen, length)
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_genet = sort(unique(sim_select_genetic$SA_Gen)),
prop_val_sign = tapply(sim_select_genetic$IndicGen,
sim_select_genetic$SA_Gen, mean),
prop_val_sign_aov = tapply(sim_select_genetic$IndicGen_aov,
sim_select_genetic$SA_Gen, mean))
#Plot
plot_genetic_poisson<-ggplot2::ggplot(data = data_prop, aes(x=val_SA_genet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 4) +
geom_line(aes(y=prop_val_sign, color="LMER", linetype = "Balanced"), size = 1.5) +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova", linetype = "Balanced"), size = 0.75) +
ylab("Proportion of simulations with\nsignificant genetic local adaption") +
xlab("Intensity of genetic local adaptation (SA)") +
ggtitle("Poisson data") +
theme_LO_sober +
labs(color = "Methods", linetype = "Design") +
scale_color_manual(values = colors) +
scale_linetype_manual(values = lines)
plot_genetic_poisson
# #Save plot
# name_plot<-paste0("Simul_powertest_genetic_poisson_",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# plot_genetic_poisson, base_height = 10/cm(1), base_width = 14/cm(1), dpi = 1200)
#################################
### Plot non-genetic
#################################
#Number of simuls for each SA non genetic value
tapply(sim$index_SA_NonGen, sim$SA_NonGen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 4 111 493 1054 1759 1912 1695 1232 853 451 267 98 52 13 5 1#Save the same number of simuls for each SA value
nb_per_SA ## [1] 100sim_select_nongenetic <- select_sample_simul(dataset_sim = sim,
type_SA = "non-genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#To check
# dim(sim)
# dim(sim_select_nongenetic)
# tapply(sim_select_nongenetic$index_SA_Gen, sim_select_nongenetic$SA_Gen, length)
# tapply(sim_select_nongenetic$index_SA_NonGen, sim_select_nongenetic$SA_NonGen, length)
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_nongenet = sort(unique(sim_select_nongenetic$SA_NonGen)),
prop_val_sign = tapply(sim_select_nongenetic$IndicNonGen,
sim_select_nongenetic$SA_NonGen, mean),
prop_val_sign_aov = tapply(sim_select_nongenetic$IndicNonGen_aov,
sim_select_nongenetic$SA_NonGen, mean))
colors<-c("LMER" = "#DCA237", "Anova" = "#459B76")
#Plot
plot_nongenetic_poisson<-ggplot2::ggplot(data = data_prop, aes(x=val_SA_nongenet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 4) +
geom_line(aes(y=prop_val_sign, color="LMER", linetype = "Balanced"), size = 1.5) +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova", linetype = "Balanced"), size = 0.75) +
ylab("Proportion of simulations with\nsignificant non-genetic local adaption") +
xlab("Intensity of non-genetic local adaptation (SA)") +
ggtitle("Poisson data") +
theme_LO_sober +
labs(color = "Methods", linetype = "Design") +
scale_color_manual(values = colors) +
scale_linetype_manual(values = lines)
plot_nongenetic_poisson
#Save plot
# name_plot<-paste0("Simul_powertest_nongenetic_poisson_",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# plot_nongenetic_poisson, base_height = 10/cm(1), base_width = 14/cm(1), dpi = 1200)
#################################
### Common plot
#################################
power_plot_poisson<-cowplot::plot_grid(plot_nongenetic_poisson,plot_genetic_poisson,
labels=c("A", "B"),
label_size = 15,
ncol =1, nrow = 2,
# hjust = 0, vjust = 1,
scale = c(1, 1))
power_plot_poisson
#Save plot
# name_plot<-paste0("Simul_powertest_poissondata_",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# power_plot_poisson, base_height = 20/cm(1), base_width = 14/cm(1), dpi = 1200)
#
# 4.3 False positive rate analyses
### No genetic and no non genetic effects
nb_simul = 200
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "poisson", rho=0, rho_ng=0,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
mean(sim$IndicGen)## [1] 0.025mean(sim$IndicGen_aov)## [1] 0.025## Local adaptation (non genetic)
mean(sim$IndicNonGen)## [1] 0.01mean(sim$IndicNonGen_aov)## [1] 0.03### Genetic effects and no non-genetic effects
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "poisson", rho=-0.1, rho_ng=0,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
tapply(sim$IndicGen, sim$SA_Gen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## 0.2500000 0.0000000 0.1333333 0.4500000 0.2564103 0.2333333 0.3846154 0.3478261
## 1 1.1 1.2 1.3
## 0.5454545 0.8333333 1.0000000 1.0000000tapply(sim$IndicGen_aov, sim$SA_Gen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## 0.2500000 0.0000000 0.1333333 0.4500000 0.2564103 0.2333333 0.3846154 0.3478261
## 1 1.1 1.2 1.3
## 0.5454545 0.8333333 1.0000000 1.0000000## Local adaptation (non genetic)
mean(sim$IndicNonGen)## [1] 0.045mean(sim$IndicNonGen_aov)## [1] 0.05### Non genetic effects and no genetic effects
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "poisson", rho=0, rho_ng=-0.1,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
mean(sim$IndicGen)## [1] 0.095mean(sim$IndicGen_aov)## [1] 0.095## Local adaptation (non genetic)
tapply(sim$IndicNonGen, sim$SA_NonGen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## 0.0000000 0.0000000 0.0500000 0.1250000 0.3673469 0.2432432 0.2916667 0.5333333
## 1 1.1 1.2
## 0.2500000 0.5000000 1.0000000tapply(sim$IndicNonGen_aov, sim$SA_NonGen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## 0.0000000 0.0000000 0.1500000 0.1875000 0.2857143 0.2702703 0.3750000 0.4000000
## 1 1.1 1.2
## 0.2500000 0.5000000 1.00000005 Binomial distribution
5.1 Simulation
nb_simul <- 10000
#Perform siluations
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "binomial", rho = -0.05,
rho_ng = -0.05,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5, ntrial = 60)
sim <- as.data.frame(t(sim))
#Add indic SA
sim <- add_indic_sign(sim)
#Add number of simul for each SA value to keep only 100 simuls per SA value
sim <- add_sim_number_SA(sim)
## Local adaptation (genetic)
tapply(sim$IndicGen[sim$index_SA_Gen<101], sim$SA_Gen[sim$index_SA_Gen<101], length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 100 100 100 100 100 100 100 100 100 100 100 37 19 1tapply(sim$IndicGen_aov[sim$index_SA_Gen<101], sim$SA_Gen[sim$index_SA_Gen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.6000000 0.6300000 0.6000000 0.5400000 0.6500000 0.6500000 0.7200000 0.7000000
## 0.9 1 1.1 1.2 1.3 1.4 1.5
## 0.6700000 0.8100000 0.8100000 0.8400000 0.8918919 0.8421053 1.0000000## Local adaptation (non genetic)
tapply(sim$IndicNonGen[sim$index_SA_NonGen<101], sim$SA_NonGen[sim$index_SA_NonGen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.0000000 0.2000000 0.1600000 0.2200000 0.2700000 0.3500000 0.3700000 0.5300000
## 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 0.5500000 0.5400000 0.6500000 0.7040816 0.6538462 0.9230769 0.8000000 1.0000000tapply(sim$IndicNonGen_aov[sim$index_SA_NonGen<101], sim$SA_NonGen[sim$index_SA_NonGen<101], mean)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.0000000 0.2000000 0.1700000 0.2100000 0.2700000 0.3500000 0.3700000 0.5400000
## 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 0.5500000 0.5400000 0.6500000 0.7040816 0.6538462 0.9230769 0.8000000 1.0000000tapply(sim$IndicGen, sim$SA_Gen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 107 496 1127 1633 1919 1672 1307 868 464 234 106 37 19 1tapply(sim$IndicNonGen, sim$SA_NonGen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 4 111 493 1054 1759 1912 1695 1232 853 451 267 98 52 13 5 15.2 Plot
####################################
###### Plot genetic
####################################
#Number of simuls for each SA value
tapply(sim$index_SA_Gen, sim$SA_Gen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 107 496 1127 1633 1919 1672 1307 868 464 234 106 37 19 1#Save the same number of simuls for each SA value
nb_per_SA=100
sim_select_genetic <- select_sample_simul(dataset_sim = sim,
type_SA = "genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#To check
# dim(sim)
# dim(sim_select_genetic)
# tapply(sim_select_genetic$index_SA_Gen, sim_select_genetic$SA_Gen, length)
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_genet = sort(unique(sim_select_genetic$SA_Gen)),
prop_val_sign = tapply(sim_select_genetic$IndicGen,
sim_select_genetic$SA_Gen, mean),
prop_val_sign_aov = tapply(sim_select_genetic$IndicGen_aov,
sim_select_genetic$SA_Gen, mean))
#Plot
plot_genetic_binomial<-ggplot2::ggplot(data = data_prop, aes(x=val_SA_genet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 4) +
geom_line(aes(y=prop_val_sign, color="LMER", linetype = "Balanced"), size = 1.5) +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova", linetype = "Balanced"), size = 0.75) +
ylab("Proportion of simulations with\nsignificant genetic local adaption") +
xlab("Intensity of genetic local adaptation (SA)") +
ggtitle("Binomial data") +
theme_LO_sober +
labs(color = "Methods", linetype = "Design") +
scale_color_manual(values = colors) +
scale_linetype_manual(values = lines)
plot_genetic_binomial
# #Save plot
# name_plot<-paste0("Simul_powertest_genetic_binomial_",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# plot_genetic_binomial, base_height = 10/cm(1), base_width = 14/cm(1), dpi = 1200)
#################################
### Plot non-genetic
#################################
#Number of simuls for each SA non genetic value
tapply(sim$index_SA_NonGen, sim$SA_NonGen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 4 111 493 1054 1759 1912 1695 1232 853 451 267 98 52 13 5 1#Save the same number of simuls for each SA value
nb_per_SA ## [1] 100sim_select_nongenetic <- select_sample_simul(dataset_sim = sim,
type_SA = "non-genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#To check
# dim(sim)
# dim(sim_select_nongenetic)
# tapply(sim_select_nongenetic$index_SA_Gen, sim_select_nongenetic$SA_Gen, length)
# tapply(sim_select_nongenetic$index_SA_NonGen, sim_select_nongenetic$SA_NonGen, length)
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_nongenet = sort(unique(sim_select_nongenetic$SA_NonGen)),
prop_val_sign = tapply(sim_select_nongenetic$IndicNonGen,
sim_select_nongenetic$SA_NonGen, mean),
prop_val_sign_aov = tapply(sim_select_nongenetic$IndicNonGen_aov,
sim_select_nongenetic$SA_NonGen, mean))
colors<-c("LMER" = "#DCA237", "Anova" = "#459B76")
#Plot
plot_nongenetic_binomial<-ggplot2::ggplot(data = data_prop, aes(x=val_SA_nongenet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 4) +
geom_line(aes(y=prop_val_sign, color="LMER", linetype = "Balanced"), size = 1.5) +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova", linetype = "Balanced"), size = 0.75) +
ylab("Proportion of simulations with\nsignificant non-genetic local adaption") +
xlab("Intensity of non-genetic local adaptation (SA)") +
ggtitle("Binomial data") +
theme_LO_sober +
labs(color = "Methods", linetype = "Design") +
scale_color_manual(values = colors) +
scale_linetype_manual(values = lines)
plot_nongenetic_binomial
#Save plot
# name_plot<-paste0("Simul_powertest_nongenetic_binomial_",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# plot_nongenetic_binomial, base_height = 10/cm(1), base_width = 14/cm(1), dpi = 1200)
#################################
### Common plot
#################################
power_plot_binomial<-cowplot::plot_grid(plot_nongenetic_binomial,plot_genetic_binomial,
labels=c("A", "B"),
label_size = 15,
ncol =1, nrow = 2,
# hjust = 0, vjust = 1,
scale = c(1, 1))
power_plot_binomial
#Save plot
# name_plot<-paste0("Simul_powertest_binomialdata_",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# power_plot_binomial, base_height = 20/cm(1), base_width = 14/cm(1), dpi = 1200)
#
# 5.3 All plots
power_plot_all<-cowplot::plot_grid(plot_nongenetic_normal,plot_nongenetic_poisson,plot_nongenetic_binomial,
plot_genetic_normal, plot_genetic_poisson,plot_genetic_binomial,
#labels=c("A", "B"),
#label_size = 15,
ncol =3, nrow = 2,
# hjust = 0, vjust = 1,
scale = c(1, 1))
power_plot_all
#Save plot
name_plot<-paste0("Simul_powertest_ALLDATA_",
nb_simul,"nbsimul_",nb_per_SA,
"nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# power_plot_all, base_height = 20/cm(1), base_width = 42/cm(1), dpi = 1200)
# 5.4 False positive rate analyses
### No genetic and no non genetic effects
nb_simul = 200
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "binomial", rho=0, rho_ng=0,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
mean(sim$IndicGen)## [1] 0.01mean(sim$IndicGen_aov)## [1] 0.01## Local adaptation (non genetic)
mean(sim$IndicNonGen)## [1] 0.015mean(sim$IndicNonGen_aov)## [1] 0.015### Genetic effects and no non-genetic effects
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "binomial",rho=-0.1, rho_ng=0,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
tapply(sim$IndicGen, sim$SA_Gen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.00000000 0.16666667 0.06666667 0.45000000 0.17948718 0.23333333 0.26923077
## 0.9 1 1.1 1.2 1.3
## 0.39130435 0.36363636 0.50000000 1.00000000 1.00000000tapply(sim$IndicGen_aov, sim$SA_Gen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8
## 0.00000000 0.16666667 0.06666667 0.45000000 0.17948718 0.23333333 0.26923077
## 0.9 1 1.1 1.2 1.3
## 0.39130435 0.36363636 0.50000000 1.00000000 1.00000000## Local adaptation (non genetic)
mean(sim$IndicNonGen)## [1] 0.015mean(sim$IndicNonGen_aov)## [1] 0.015###
### Non genetic effects and no genetic effects
sim <- sapply(1:nb_simul, simul_fitnessdata, distrib = "binomial",rho=0, rho_ng=-0.1,
npop_per_fruit = 5, nfruit = 3, nhab = 3, nrep = 10, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
## Local adaptation (genetic)
mean(sim$IndicGen)## [1] 0.07mean(sim$IndicGen_aov)## [1] 0.07## Local adaptation (non genetic)
tapply(sim$IndicNonGen, sim$SA_NonGen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## 0.0000000 0.2500000 0.1000000 0.1875000 0.3061224 0.3783784 0.5000000 0.6000000
## 1 1.1 1.2
## 0.2500000 0.5000000 1.0000000tapply(sim$IndicNonGen_aov, sim$SA_NonGen, mean)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## 0.0000000 0.2500000 0.2000000 0.1875000 0.3061224 0.3783784 0.5000000 0.6000000
## 1 1.1 1.2
## 0.2500000 0.5000000 1.00000006 UNBALANCED DESIGN
#Dataset
#unbalanced dataset
data_PERF <- read.table(file="/Users/laure/Documents/RESEARCH/z_PhD/POPNAT/OLD/DATA_COMPLET/Clean/DATACOMPLET_PERF_POPNAT2018.csv",sep=";",header=T)
data_PERF$Original_environment<-as.factor(data_PERF$Original_environment)
data_PERF <- data_PERF[data_PERF$Original_environment!="WT3",]
data_PERF <- data_PERF[data_PERF$Test_environment!="GF",]
data_PERF <- data_PERF[data_PERF$Test_environment!="Grape",] #mvrnom doesn't work with unbalanced matrix
data_PERF <- droplevels(data_PERF)7 Normal distribution
7.1 Simulation
#Simul data
nb_simul <- 10000
sim <- sapply(1:nb_simul, simul_fitnessdata_unbalanced, distrib = "normal",
unbalanced_dataset = data_PERF, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5, rho = -0.05, rho_ng = -0.05,
sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
sim <- add_sim_number_SA(sim)7.2 Plot
####################################
###### Plot genetic
####################################
#Number of simuls for each SA value
tapply(sim$index_SA_Gen, sim$SA_Gen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 107 496 1127 1633 1919 1672 1307 868 464 234 106 37 19 1#Save the same number of simuls for each SA value
nb_per_SA=100
sim_select_genetic <- select_sample_simul(dataset_sim = sim,
type_SA = "genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_genet = sort(unique(sim_select_genetic$SA_Gen)),
prop_val_sign = tapply(sim_select_genetic$IndicGen,
sim_select_genetic$SA_Gen, mean),
prop_val_sign_aov = tapply(sim_select_genetic$IndicGen_aov,
sim_select_genetic$SA_Gen, mean))
#Plot
plot_genetic_normal_unbalanced<-ggplot2::ggplot(data = data_prop, aes(x=val_SA_genet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 2.5) +
geom_line(aes(y=prop_val_sign, color="LMER"), size = 1, linetype = "dotted") +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova"), size = 1, linetype = "dotted") +
ylab("Proportion of simulations with\nsignificant genetic local adaption") +
xlab("Intensity of genetic local adaptation (SA)") +
ggtitle("Normal data (unbalanced experimental design)") +
theme_LO_sober +
labs(color = "Methods", linetype = "Design") +
scale_color_manual(values = colors) +
scale_linetype_manual(values = lines)
plot_genetic_normal_unbalanced
# Balanced + unbabacaled
plot_genetic_normal_ALL <- plot_genetic_normal +
geom_point(data = data_prop,
aes(y=prop_val_sign, x=val_SA_genet, color="LMER"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign, x=val_SA_genet,
color="LMER", linetype = "Unbalanced"), size = 1) +
geom_point(data = data_prop,
aes(y=prop_val_sign_aov, x=val_SA_genet,
color="Anova"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign_aov, color="Anova", x=val_SA_genet,
linetype = "Unbalanced"), size = 1) +
labs(color = "Methods", linetype = "Design") +
scale_linetype_manual(values = lines)
plot_genetic_normal_ALL
#################################
### Plot non-genetic
#################################
#Number of simuls for each SA non genetic value
tapply(sim$index_SA_NonGen, sim$SA_NonGen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 4 111 493 1054 1759 1912 1695 1232 853 451 267 98 52 13 5 1#Save the same number of simuls for each SA value
nb_per_SA ## [1] 100sim_select_nongenetic <- select_sample_simul(dataset_sim = sim,
type_SA = "non-genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_nongenet = sort(unique(sim_select_nongenetic$SA_NonGen)),
prop_val_sign = tapply(sim_select_nongenetic$IndicNonGen,
sim_select_nongenetic$SA_NonGen, mean),
prop_val_sign_aov = tapply(sim_select_nongenetic$IndicNonGen_aov,
sim_select_nongenetic$SA_NonGen, mean))
#Plot
plot_nongenetic_normal_unbalanced<-ggplot2::ggplot(data = data_prop,
aes(x=val_SA_nongenet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 2.5) +
geom_line(aes(y=prop_val_sign, color="LMER"), size = 1, linetype = "dotted") +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova"), size = 1, linetype = "dotted") +
ylab("Proportion of simulations with\nsignificant non-genetic local adaption") +
xlab("Intensity of non-genetic local adaptation (SA)") +
theme_LO_sober +
ggtitle("Normal data (unbalanced experimental design)") +
labs(color = "Methods") +
scale_color_manual(values = colors)
plot_nongenetic_normal_unbalanced
plot_nongenetic_normal_ALL <- plot_genetic_normal +
geom_point(data = data_prop,
aes(y=prop_val_sign, x=val_SA_nongenet, color="LMER"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign, x=val_SA_nongenet,
color="LMER", linetype = "Unbalanced"), size = 1) +
geom_point(data = data_prop,
aes(y=prop_val_sign_aov, x=val_SA_nongenet,
color="Anova"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign_aov, color="Anova", x=val_SA_nongenet,
linetype = "Unbalanced"), size = 1) +
labs(color = "Methods", linetype = "Design") +
scale_linetype_manual(values = lines)
plot_nongenetic_normal_ALL
8 Poisson distribution
8.1 Simulation
#Simul data
#nb_simul <- 1000
sim <- sapply(1:nb_simul, simul_fitnessdata_unbalanced, distrib = "poisson",
unbalanced_dataset = data_PERF, sdpop = 0.5,
sdfruithab = 1, sdfruithab_ng = 1, rho = -0.05, rho_ng = -0.05, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
sim <- add_sim_number_SA(sim)8.2 Plot
####################################
###### Plot genetic
####################################
#Number of simuls for each SA value
tapply(sim$index_SA_Gen, sim$SA_Gen, length)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2 2.1
## 3 10 50 121 253 404 558 709 819 939 955 915 814 777 652 531 465 313 221 171
## 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3.1
## 119 63 60 36 15 12 11 3 1#Save the same number of simuls for each SA value
nb_per_SA=100
sim_select_genetic <- select_sample_simul(dataset_sim = sim,
type_SA = "genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_genet = sort(unique(sim_select_genetic$SA_Gen)),
prop_val_sign = tapply(sim_select_genetic$IndicGen,
sim_select_genetic$SA_Gen, mean),
prop_val_sign_aov = tapply(sim_select_genetic$IndicGen_aov,
sim_select_genetic$SA_Gen, mean))
#Plot
plot_genetic_poisson_unbalanced<-ggplot2::ggplot(data = data_prop, aes(x=val_SA_genet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 2.5) +
geom_line(aes(y=prop_val_sign, color="LMER"), size = 1, linetype = "dotted") +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova"), size = 1, linetype = "dotted") +
ylab("Proportion of simulations with\nsignificant genetic local adaption") +
xlab("Intensity of genetic local adaptation (SA)") +
ggtitle("Poisson data (unbalanced experimental design)") +
theme_LO_sober +
labs(color = "Methods") +
scale_color_manual(values = colors)
plot_genetic_poisson_unbalanced
plot_genetic_poisson_ALL <- plot_genetic_poisson +
geom_point(data = data_prop,
aes(y=prop_val_sign, x=val_SA_genet, color="LMER"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign, x=val_SA_genet,
color="LMER", linetype = "Unbalanced"), size = 1) +
geom_point(data = data_prop,
aes(y=prop_val_sign_aov, x=val_SA_genet,
color="Anova"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign_aov, color="Anova", x=val_SA_genet,
linetype = "Unbalanced"), size = 1) +
labs(color = "Methods", linetype = "Design") +
scale_linetype_manual(values = lines)
plot_genetic_poisson_ALL
#################################
### Plot non-genetic
#################################
#Number of simuls for each SA non genetic value
tapply(sim$index_SA_NonGen, sim$SA_NonGen, length)## 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2 2.1
## 2 11 49 129 240 379 522 727 885 969 951 905 845 743 663 496 405 325 224 183
## 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3 3.3
## 127 78 53 35 26 15 5 3 4 1#Save the same number of simuls for each SA value
nb_per_SA ## [1] 100sim_select_nongenetic <- select_sample_simul(dataset_sim = sim,
type_SA = "non-genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_nongenet = sort(unique(sim_select_nongenetic$SA_NonGen)),
prop_val_sign = tapply(sim_select_nongenetic$IndicNonGen,
sim_select_nongenetic$SA_NonGen, mean),
prop_val_sign_aov = tapply(sim_select_nongenetic$IndicNonGen_aov,
sim_select_nongenetic$SA_NonGen, mean))
colors<-c("LMER" = "#DCA237", "Anova" = "#459B76")
#Plot
plot_nongenetic_poisson_unbalanced<-ggplot2::ggplot(data = data_prop,
aes(x=val_SA_nongenet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 2.5) +
geom_line(aes(y=prop_val_sign, color="LMER"), size = 1, linetype = "dotted") +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova"), size = 1, linetype = "dotted") +
ylab("Proportion of simulations with\nsignificant non-genetic local adaption") +
xlab("Intensity of non-genetic local adaptation (SA)") +
theme_LO_sober +
ggtitle("Poisson data (unbalanced experimental design)") +
labs(color = "Methods") +
scale_color_manual(values = colors)
plot_nongenetic_poisson_unbalanced
plot_nongenetic_poisson_ALL <- plot_nongenetic_poisson +
geom_point(data = data_prop,
aes(y=prop_val_sign, x=val_SA_nongenet, color="LMER"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign, x=val_SA_nongenet,
color="LMER", linetype = "Unbalanced"), size = 1) +
geom_point(data = data_prop,
aes(y=prop_val_sign_aov, x=val_SA_nongenet,
color="Anova"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign_aov, color="Anova", x=val_SA_nongenet,
linetype = "Unbalanced"), size = 1) +
labs(color = "Methods", linetype = "Design") +
scale_linetype_manual(values = lines)
plot_nongenetic_poisson_ALL
9 Binomial distribution
9.1 Simulation
#Simul data
#nb_simul <- 1000
sim <- sapply(1:nb_simul, simul_fitnessdata_unbalanced, distrib = "binomial",
unbalanced_dataset = data_PERF, sdpop = 0.5,
sdfruithab = 0.5, sdfruithab_ng = 0.5,
rho = -0.05, rho_ng = -0.05, sigma = 0.5)
sim <- as.data.frame(t(sim))
sim <- add_indic_sign(sim)
sim <- add_sim_number_SA(sim)9.2 Plot
####################################
###### Plot genetic
####################################
#Number of simuls for each SA value
tapply(sim$index_SA_Gen, sim$SA_Gen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5
## 10 107 496 1127 1633 1919 1672 1307 868 464 234 106 37 19 1#Save the same number of simuls for each SA value
nb_per_SA=100
sim_select_genetic <- select_sample_simul(dataset_sim = sim,
type_SA = "genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_genet = sort(unique(sim_select_genetic$SA_Gen)),
prop_val_sign = tapply(sim_select_genetic$IndicGen,
sim_select_genetic$SA_Gen, mean),
prop_val_sign_aov = tapply(sim_select_genetic$IndicGen_aov,
sim_select_genetic$SA_Gen, mean))
#Plot
plot_genetic_Binomial_unbalanced<-ggplot2::ggplot(data = data_prop,
aes(x=val_SA_genet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 2.5) +
geom_line(aes(y=prop_val_sign, color="LMER"), size = 1, linetype = "dotted") +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova"), size = 1, linetype = "dotted") +
ylab("Proportion of simulations with\nsignificant genetic local adaption") +
xlab("Intensity of genetic local adaptation (SA)") +
ggtitle("Binomial data (unbalanced experimental design)") +
theme_LO_sober +
labs(color = "Methods") +
scale_color_manual(values = colors)
plot_genetic_Binomial_unbalanced
plot_genetic_binomial_ALL <- plot_genetic_binomial +
geom_point(data = data_prop,
aes(y=prop_val_sign, x=val_SA_genet, color="LMER"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign, x=val_SA_genet,
color="LMER", linetype = "Unbalanced"), size = 1) +
geom_point(data = data_prop,
aes(y=prop_val_sign_aov, x=val_SA_genet,
color="Anova"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign_aov, color="Anova", x=val_SA_genet,
linetype = "Unbalanced"), size = 1) +
labs(color = "Methods", linetype = "Design") +
scale_linetype_manual(values = lines)
plot_genetic_binomial_ALL
#################################
### Plot non-genetic
#################################
#Number of simuls for each SA non genetic value
tapply(sim$index_SA_NonGen, sim$SA_NonGen, length)## 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 1.7
## 4 111 493 1054 1759 1912 1695 1232 853 451 267 98 52 13 5 1#Save the same number of simuls for each SA value
nb_per_SA ## [1] 100sim_select_nongenetic <- select_sample_simul(dataset_sim = sim,
type_SA = "non-genetic",
nb_simul_per_SA = nb_per_SA,
minimum_SA = "0")
#Compute proportion of significant simul for each SA values
data_prop <- data.frame(val_SA_nongenet = sort(unique(sim_select_nongenetic$SA_NonGen)),
prop_val_sign = tapply(sim_select_nongenetic$IndicNonGen,
sim_select_nongenetic$SA_NonGen, mean),
prop_val_sign_aov = tapply(sim_select_nongenetic$IndicNonGen_aov,
sim_select_nongenetic$SA_NonGen, mean))
colors<-c("LMER" = "#DCA237", "Anova" = "#459B76")
#Plot
plot_nongenetic_Binomial_unbalanced<-ggplot2::ggplot(data = data_prop,
aes(x=val_SA_nongenet)) +
geom_point(aes(y=prop_val_sign, color="LMER"), size = 2.5) +
geom_line(aes(y=prop_val_sign, color="LMER"), size = 1, linetype = "dotted") +
geom_point(aes(y=prop_val_sign_aov, color="Anova"), size = 2.5) +
geom_line(aes(y=prop_val_sign_aov, color="Anova"), size = 1, linetype = "dotted") +
ylab("Proportion of simulations with\nsignificant non-genetic local adaption") +
xlab("Intensity of non-genetic local adaptation (SA)") +
theme_LO_sober +
ggtitle("Binomial data (unbalanced experimental design)") +
labs(color = "Methods") +
scale_color_manual(values = colors)
plot_nongenetic_Binomial_unbalanced
plot_nongenetic_binomial_ALL <- plot_nongenetic_binomial +
geom_point(data = data_prop,
aes(y=prop_val_sign, x=val_SA_nongenet, color="LMER"),
size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign, x=val_SA_nongenet,
color="LMER", linetype = "Unbalanced"), size = 1) +
geom_point(data = data_prop,
aes(y=prop_val_sign_aov, x=val_SA_nongenet,
color="Anova"), size = 2.5) +
geom_line(data = data_prop,
aes(y=prop_val_sign_aov, color="Anova", x=val_SA_nongenet,
linetype = "Unbalanced"), size = 1) +
labs(color = "Methods", linetype = "Design") +
scale_linetype_manual(values = lines)
plot_nongenetic_binomial_ALL
10 All plots unbalanced
power_plot_all_unbalanced<-cowplot::plot_grid(plot_nongenetic_normal_unbalanced,
plot_nongenetic_poisson_unbalanced,
plot_nongenetic_Binomial_unbalanced,
plot_genetic_normal_unbalanced,
plot_genetic_poisson_unbalanced,
plot_genetic_Binomial_unbalanced,
#labels=c("A", "B"),
#label_size = 15,
ncol =3, nrow = 2,
# hjust = 0, vjust = 1,
scale = c(1, 1))
power_plot_all_unbalanced
#Save plot
name_plot<-paste0("Simul_powertest_ALLDATA_UNBALANCED",
nb_simul,"nbsimul_",nb_per_SA,
"nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# power_plot_all_unbalanced, base_height = 20/cm(1),
# base_width = 42/cm(1), dpi = 1200)
# 11 ALLPLOT
POWER_ALL<-cowplot::plot_grid(plot_nongenetic_normal,
plot_nongenetic_poisson,
plot_nongenetic_binomial,
plot_genetic_normal,
plot_genetic_poisson,
plot_genetic_binomial,
plot_nongenetic_normal_unbalanced,
plot_nongenetic_poisson_unbalanced,
plot_nongenetic_Binomial_unbalanced,
plot_genetic_normal_unbalanced,
plot_genetic_poisson_unbalanced,
plot_genetic_Binomial_unbalanced,
#labels=c("A", "B"),
#label_size = 15,
ncol =3, nrow = 4,
# hjust = 0, vjust = 1,
scale = c(1, 1))
POWER_ALL
#Save plot
# name_plot<-paste0("Simul_powertest_ALLDATA_BALANCED_and_UNBALANCED",
# nb_simul,"nbsimul_",nb_per_SA,
# "nbsimulperSA.pdf")
# cowplot::save_plot(file =here::here("figures",name_plot ),
# POWER_ALL, base_height = 40/cm(1),
# base_width = 42/cm(1), dpi = 1200)
#
PLOT_ALL<-cowplot::plot_grid(plot_nongenetic_normal_ALL,
plot_nongenetic_poisson_ALL,
plot_nongenetic_binomial_ALL,
plot_genetic_normal_ALL,
plot_genetic_poisson_ALL,
plot_genetic_binomial_ALL,
#labels=c("A", "B"),
#label_size = 15,
ncol =3, nrow = 2,
# hjust = 0, vjust = 1,
scale = c(1, 1))
PLOT_ALL
#Save plot
name_plot<-paste0("Simul_powertest_",
nb_simul,"nbsimul_",nb_per_SA,
"nbsimulperSA.pdf")
cowplot::save_plot(file = here::here("figures",name_plot ),
PLOT_ALL, base_height = 20/cm(1),
base_width = 42/cm(1), dpi = 1200)